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Abstract 

A binary embedded-atom method (EAM) potential is optimized for Cu on Ag(ll 1) by fitting to ab initio 
data. The fitting database consists of DFT calculations of Cu monomers and dimers on Ag(l 11), specifically 
their relative energies, adatom heights, and dimer separations. We start from the Mishin Cu-Ag EAM 
potential and first modify the Cu-Ag pair potential to match the FCC/HCP site energy difference then include 
Cu-Cu pair potential optimization for the entire database. The optimized EAM potential reproduce DFT 
monomer and dimer relative energies and geometries correctly. In trimer calculations, the potential produces 
the DFT relative energy between FCC and HCP trimers, though a different ground state is predicted. We use 
the optimized potential to calculate diffusion barriers for Cu monomers, dimers, and trimers. The predicted 
monomer barrier is the same as DFT, while experimental barriers for monomers and dimers are both lower 
than predicted here. We attribute the difference with experiment to the overestimation of surface adsorption 
energies by DFT and a simple correction is presented. Our results show that the optimized Cu-Ag EAM can 
be applied in the study of larger Cu islands on Ag(l 1 1). 

PACS numbers: 68.35.Fx,68.35.bd,68.35.Gy,71.15.Pd 
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I. INTRODUCTION 



Knowledge of the surface diffusion dynamics for small atom clusters is critical to understanding 
heteroepitaxial thin film growth. While numerous experiment o 1 ' 2 ' 3 ' 4 and computer simulations^^ 
have studied homogeneous systems, less is known about lattice mismatched heterogeneous 
systems 9 *^2iii and their interesting diffusion kinetics. In this study, we consider Cu on Ag(l 1 1),— 
a system with a lattice mismatch of 12%. 13 The lattice mismatch induces strain in both the island 
and substrate and has been predicted to promote rapid diffusion. 14 

To accurately compute the energetics of surface island systems, first principle density- 
functional theory (DFT) calculations are preferred to empirical potentials. However, DFT methods 
are too computationally intensive to efficiently search the phase space of each island and accu- 
rate classical potentials are needed to characterize island diffusional dynamics. The embedded 
atom method 15 (EAM) is well suited for metallic systems combining pair interactions with an 
atomic embedding energy term dependent on the local "electron density." Table [I] shows that other 
EAM potentials were unable to reproduce DFT calculated Cu island energies and geometries on 
Ag(l 1 1), motivating the search for a new potential. 

We optimize a new EAM potential for Cu on Ag(lll) using monomer and dimer data. Sec- 
tion |n] explains the DFT and EAM calculation parameters in detail. Section [HI] presents the pro- 
cedure for the potential optimization. The energetics and diffusion results from the new EAM for 
monomers, dimers, and trimers are reported in Section [IV] We justify the new potential for the 
study of small Cu islands on Ag(l 1 1) surface by comparing the calculated results to experimental 
and DFT values in Section |Vl 

II. COMPUTATIONAL DETAILS 

The density-functional theory calculations are performed with vasp,^*^ a density-functional 
code using a plane-wave basis and ultrasoft Vanderbilt-type pseudopotentials.^2 The local- 
density approximation as parametrized by Perdew and Zunger 20 and a plane-wave kinetic-energy 
cutoff of 200 eV ensures accurate treatment of the Cu and Ag potential. We treat the s and d states 
as valence, corresponding to an Ar and Kr core atomic reference configuration for Cu and Ag, 
respectively. The (111) surface slab calculations used a 3 x 3 geometry with 6(111) planes of Ag 
and 6 (111) planes of vacuum; the fc-point meshes for the surface slab calculations are 8 x 8 x 1, 
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with a Methfessel-Paxton smearing of 0.25 eV. 

EAM energy values were computed with the lammps molecular dynamics package, 21 The 
monomer and dimer results in Table U are obtained using a periodic 3x3 cell of 6 (111) planes. 
The trimers are calculated with 4x4 periodic cells. Results presented in section IV are from 6x6 
periodic cells, where our potential predicts a finite-size effect of less than 5meV compared to the 
3x3 cell. Transition energy barriers are determined with nudged elastic band 22 calculations after 
initial and final states have been found through molecular dynamics or the dimer search method.— 
Attempt frequency prefactors are computed with the Vineyard formula, 24 taking the ratio between 
the product of harmonic vibrational frequencies at the initial state and the saddle point. 

DDL OPTIMIZATION PROCEDURE 

In EAM, the total energy of the system is given by 



where 0,7(^7) is the pair potential interaction between atoms i and j separated by a distance of r (J - 
and Fi(p im ) is embedding energy of atom i in the superposition of atomic electron densities Pj(jij). 
The Mishin CuAg binary EAM potential 2 ^! is described by seven functions: 0euCu(/)> 0AgAg(/)> 
0CuA g (/X Pcu(r), PA g (r), F Cu (p), and F Ag (p). The Mishin EAM embedding energy functions and 
electron density functions are not changed in our optimization. Only the Cu-Ag and Cu-Cu pair 
potentials are modified to fit our DFT optimization database. We forgo modification of the Ag-Ag 
potential because the distance between relaxed EAM Ag(l 1 1) planes are within 3% of the relaxed 
DFT Ag surface. 

Cu monomers and dimers are building blocks for larger islands, making them ideal choices for 
the optimization database. The Ag(lll) surface is divided into FCC and HCP sites, depending 
on the atomic configuration continuing from the top two layers of Ag. For monomers, the single 
Cu atom rests at either a FCC or HCP site. For dimers, four different configuration of the Cu 
pair can be formed, FCC-FCC (FF), HCP-HCP (HH), and two types of FCC-HCP (FH short and 
FHi ong ). The two FH dimers (c.f. Fig. [5]) are differentiated by their neighboring Ag atoms, the two 
triangles of Ag neighbors can share a side (FH short ), or share a corner (FHi ong ). The optimization 
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database consists of the relative DFT energies between FCC and HCP monomers and all four 
dimers, geometric information on the heights of monomers and dimers above the Ag surface, and 
the Cu-Cu separation length. We minimize the total root-mean-square (rms) error of the energy 
differences and balance that with the total rms error of the heights and lengths. 




Height from surface [A] 



FIG. 1: Force on a Cu atom above an FCC site or an HCP site on the unrelaxed Ag(l 11) surface calculated 
by DFT, Mishin EAM, and the new optimized EAM. The ultrasoft pseudopotential DFT force displays a 
stronger and wider interaction between the Cu atom and the Ag surface than Mishin EAM. The average 
force difference (FCC and HCP) between DFT and Mishin EAM is used to produce the force-matched 
0CuAg( r ) i n Fig. 12 The final optimized EAM potential deviates between 3A-5A with a maximum deviation 
of 0. 15eV/A with respect to the DFT calculated forces. 

In Fig.[T]the DFT force of a Cu atom evaporating from a perfect (unrelaxed) Ag(l 11) surface is 
plotted, and is used in addition to the database. Starting from a height of ^[111], the force on the 
Cu atom is computed in steps of ^[11 1] for 13 points, then in steps of ^ [111] until [111], where 
the DFT force dropped to zero. The Cu atom is directly above FCC and HCP sites to cancel forces 
in the (111) plane. The difference in force between FCC and HCP is less than 0.04eV/A for both 
DFT (max deviation at ^[111]) and EAM (max deviation at |[1 1 1]). Fig. [Qshows that DFT has a 
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stronger binding of Cu to the Ag surface than the Mishin EAM. Also plotted in Fig. Q] is the force 
calculated with our optimized EAM, which captures the deeper and wider well of DFT forces. 
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FIG. 2: The Cu-Ag pair potential at sequential steps in the optimization process. The integrated force differ- 
ence from Fig. [TJ adds onto the original Mishin ^cuAgW to produce the force-matched (fcuAgW- Imposing 
a smooth cut-off at 6A and adding a spline (knots at circles) to the force-matched potential produces the 
optimized 0cuAg(/)- The spline knots are at 3A, 3.5A, 4A, 4.5A, 5A, and 6A. 

Fig.[2]shows the Cu-Ag pair potential extracted from force-matching to the DFT data. Starting 
from the highest point and moving towards the surface, the Cu atom feels the sum of forces from 
different shells of Ag atoms within a 7.2A cut-off radius. We chose this cutoff radius because the 
DFT forces goes to zero at z=7A ([111]). For each Zk, there are m = l...m k shells, in which there 
are n m Ag atoms at distance r^ m with directional component Q, m = . The z-component of the 



force at height zu is 



F z (z k ) 



-dEu 



-dE tnt dr 



dz .f- 1 dr dz 

all atoms 

m k 



(1) 



m=l 



where (p' Cu Ag( r k,m) is the radial derivative of the pair potential. We build the function (f>' CuAg (r) as 
a cubic spline with knot points r]™ ot = mm{r k m } for each k. Starting from largest Zk to smallest, 



Eqn. Q] is solved for ^cuAg^^t" ') usm §' as needed, interpolated values of (p'cuhtM^ f° r r > r , 



,knot 



CuAg^ 



The equations are solved successively until <^ U A g ls self-consistent. A final self-consistency loop 
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over all FCC and HCP forces is performed, alternating in sequence, obtaining (f>' CuAg (r) for r in the 
range from 2.04A to 7.2A. Integrating 4>' CuKg (r) generates a quartic spline, the 0cuAg( r ) plotted in 
Fig. [21 This force-matched </>cuA g (/) possess a deeper and wider energy well, capturing the stronger 
Cu-Ag interaction from DFT. For r values smaller than 2.04A, we linearly extrapolate 0cuA g ( r )- 

The force-matched 0cuAg( r ) is refined by fitting to the monomer and dimer database. The force- 
matched 0cuA g (/) has inaccurate energies for monomers and dimers, with the HCP site 4meV be- 
low the FCC site. Modifying 0cuCu(/) does not affect monomer energies, and we find that the dimer 
energy difference between FF and HH changes by less than 5meV with the 0cuCu( r ) modifications 
we present later. We optimize the Cu-Ag pair potential with respect to monomer and homogeneous 
dimer site energy differences as the next step. We reduce the interaction range to 6 A by shifting 
the potential up by 0(5. 75 A) and using quartic splines from 5 A to 6A. The quartic splines have 
two equal spaced knots within the interval and matches the value, first and second derivatives at 
5A, and at 6A goes to zero with zero slope and zero second derivative. To differentiate between 
FCC and HCP sites, we modify the Cu-Ag interaction in the range of the second and third nearest 

oo o 

neighbors for a Cu atom on the surface by adding a cubic spline, with knots at 3.5A, 4A, 4.5A, 
5A, and fixed end points at 3A and 6A. We generate (2x5+l) 4 = 14641 possible potentials with 
different values at each knot point in steps of +20meV; optimization continues using narrower 
ranges down to +lmeV. For each sweep, we select potentials with the smallest rms monomer and 
homogeneous dimer energy errors while also selecting for quantitatively low rms Cu height errors 
and potentials without multiple minimums. In Fig.[2l this optimized 4>cuKg{r) exhibits a wider well 
than the force-matched pair potential. 

Fig. |3] shows that the optimized 0cuCu(/) gives shorter and weaker bonding between Cu atoms 
on the Ag surface than in the Mishin EAM bulk Cu. We scale the original Mishin 0cuCu( r ) in 1% 
steps from 80% to 120%, and translate in 0.01A steps from -0.15A to 0.15A; potentials with Cu 
lattice parameter outside of ±5% of the bulk value are removed. A 82% scaling and a -0.13 A 
translation reproduces all relative energy differences with a final 0.5meV range optimization of 
the 0cuAg(f)- We found during optimization that although it was possible to obtain 0.0 12A rms 
dimer separation error or 0.5meV rms energy error, these two errors grew opposite one another. 
We selected for lower energy error at the expense of geometric agreement. 

Thus, our optimization procedure with respect to DFT Cu monomers and dimers follows: (1) 
modify the (pcuA g ( r ) to match the force of an evaporating Cu atom from Ag(lll) calculated in 
DFT. (2) Reduce the interaction length of the force-matched 0cuAg(r) an d add splines to reproduce 
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radial distance r [A] 



FIG. 3: The Cu-Cu pair potential from Mishin EAM and after optimization. The original Mishin 0cuCu( r ) 
is transformed with different scalings and translations simultaneously with the optimization of the ^cuAgW- 
The optimized (fcuCuW is obtained from a scaling of 82% and a translation of -0. 13 A. This indicates shorter 
and weaker bonding between Cu atoms on the Ag surface than in Cu bulk. 

the DFT FCC/HCP site energy difference. (3) Scale and translate the 0cuCu( r ) to produce better 
energetic agreement. (4) Polish the optimized 0c U A g </) potential with the optimized 0cuCu(/). The 
final 0cuA g <>) is plotted in Fig. |2j and (fcuCuO) in Fig. |3j 

In Table HI comparison with other EAM potentials show that the optimized EAM from this 
work has better agreement to DFT calculations. Among the earlier potentials, the Mishin EAM^^ 
comes closest to the DFT energies when compared to the Foiles-Baskes-Daw (FBD) EAM 15 and 
the Voter-Chen (VC) potentials. The FBD and VC EAM potentials did not indicate any 

site energy difference between FCC and HCP. None of the earlier potentials were able to capture 
the correct DFT energy difference between FH short and FH lon? dimers. Both FBD and Mishin 
EAM came within lOmeV of the experimental monomer and dimer diffusion barriers, while the 
optimized EAM potential overestimates the experimental barriers but correctly predicts the DFT 
monomer diffusion barrier. The increased Cu-Ag interaction of the optimized potential were able 
to pull the Cu atoms closer to the Ag surface, reducing the rms height error over other potentials. 
While the earlier potentials do not come close to the DFT trimer energies, the optimized EAM 
is able to capture the correct energies for AE(Hj? ot ,Fj? ot ) and AE(H^ on ,F^ on ). The deviation for the 
trimer ground state will be discussed in Section |Vl 
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TABLE I: Results for the optimized potential and other EAM potentials compared to DFT and experi- 
mental values. Three EAM potentials, Foiles-Baskes-Daw, Voter-Chen, and Mishin were compared to the 
optimized potential from this work with respect to monomer and dimer energies, diffusion barriers and ge- 
ometries. Trimer energies for DFT and our optimized potential are also presented. The symbols indicate 
less than 20% error or less than 0. 1 A error (V), between 20% and 50% error ( — ), and 50% error or greater 
than 0.1 A error (x). AE(A,B) indicates the energy difference Ea-Eb, and E a (A— >B) indicates the activa- 
tion energy to transition from state A to state B. The ab initio results, from ultrasoft pseudopotential DFT, 
are used as the standard when available. The optimized EAM potential in this work is fit to DFT values 
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IV. RESULTS 



Fig. H] shows the geometries, relative energies, and transition barriers of Cu monomers calcu- 
lated with the optimized EAM. The 14meV energy difference between the FCC and HCP site also 
represents the difference between the transition barriers. The two transitions possible are the F— >H 
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FIG. 4: Cu monomers FCC (F) and HCP (H), relative energy differences, transition energies, and attempt 
frequency prefactors. The FCC site is the ground state and can diffuse to one of three equivalent HCP sites 
nearby; similarly the HCP site diffuses to one of three FCC sites. The rate-limiting step in the diffusion 
process is the F— >H transition with an energy barrier of 93meV. 

with a 93meV barrier and the 79meV barrier H— >F transition. The F— >H barrier is higher than the 
experimental value of 65±9meV, 12 but matches our DFT calculations for the bridging site with an 
energy difference of 96meV. DFT is known to overestimate surface adsorption, 28 and we discuss 
strategies to compensate in Section |V] The agreement with DFT is a confirmation of our potential 
since the bridging site energy is not part of the optimization database. 
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FIG. 5: Cu dimers FCC-FCC (FF), FCC-HCP neighboring (FH short ), FCC-HCP non-neighboring (FHi ong ), 
and HCP-HCP (HH), relative energy differences, transition energies, and attempt frequency prefactors. The 
FF dimer is the ground state and diffusion to the HH site is achieved through one of two FH meta-stable 
sites. The pathway through the FH s h ort dominates the FF<->HH diffusion, giving an overall rate-limiting 
barrier of 88meV (FF— >FH s h rt— >HH, 80meV - 2meV + lOmeV). 
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Fig. [5] shows the geometries, relative energies, and transition barriers of Cu dimers calculated 
with the optimized EAM. The FF dimer is the ground state and the HH dimer is 26meV higher 
in energy, about twice the monomer energy difference. The FH short and FH long dimers are two 
metastable configurations which are 78meV and 130meV higher in energy than FF, respectively. 
Dimer diffusion is more complex than that for monomers, with two intermediate states between 
FF and HH plus dimer rotation. With low barrier (<lmeV) transitions out of the FHi ong state, 
the diffusion pathway through FHi ong has a 130meV barrier for FF— >FHi ong — >HH, and 103meV 
barrier for HH— >FH long — >FF. The other diffusion pathway is more complicated, since an FH short 
dimer is more likely to transition to FF (2meV barrier) than to HH (lOmeV barrier). This results 
in a 88meV barrier for FF— >FH short — >HH (88meV - 2meV + lOmeV), and 62meV barrier for 
HH— >FH sh01t — >FF. The calculated barriers are higher than the experimental barrier of 73meV, 12 
again consistent with overestimated adsorption energies by DFT 




FIG. 6: Cu trimers FFF centered on Ag (Fj? ot ), HHH centered on Ag (H^ ot ), HHH centered on a hole (H^ on ), 
and FFF centered on a hole (F^ on ), relative energy differences, transition energies, and attempt frequency 
prefactors. Fj? ot is the ground state but is only 9meV lower in energy than the F^ on state. A rotation transition 
exists between F^ ot and HjL, while the system can diffuse by overcoming the higher 289meV and 287meV 
energy barriers for F 3 ot — >H„ on and Hj? ot — >F„ on respectively. 

Fig. |6] shows the geometries, relative energies, and transition barriers of Cu trimers calculated 
with the optimized EAM. There are two different configurations for each of the FCC and HCP 
trimers due to the geometry of the (111) surface. The trimer triangles can either be centered around 
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a surface Ag atom permitting rotation, F 3 0t (groundstate) and H 3 ot (+45meV), or not, F 3 on (+9meV) 
and H 3 on (+42meV). The relative energy difference between F 3 and H 3 trimers is approximately 
two to three times the monomer energy difference. The rotation transition F 3 ot — >H 3 ot has a 149meV 
barrier and a 104meV barrier for the reverse. The non trimers do not rotate, and transition to rot 
trimers on the opposite sites. These transition barriers are higher than the rotation barriers, at 
~290meV from F 3 on — >H 3 ot and F 3 ot — >H 3 on , and ~250meV for the reverse. We expect, as with 
monomers and dimers, that the EAM overestimates the trimer transition barriers. 

We construct analytical expressions for the diffusion constants of monomers, dimers, and 
trimers using the calculated transition barriers and attempt frequencies. The rate of jumping from 
a F site to a particular H site is r FH = vf^h exp (-E a (F — » H)/k^T) where E a (F — » H) and vf^h are 
the energy barrier and the attempt frequency for the F to H transition. Then a monomer moving 
from one F to a new F site through an H site at temperature T occurs with mean wait time of 

1 monomer r, 

including the three equivalent hopping sites for each monomer transition, and with a correlation 
factor of | for monomer transitions to the original site. The Einstein diffusion relation, D = 
\ a nn T ~ l where a nn = 2.89A is the nearest-neighbor distance between Ag atoms gives the monomer 
diffusion constant 

-Dmonomer = [( r Fh) 1 + ( r Hf) 1 ] 

For both the dimer and trimer case, the diffusion system becomes complex and we use the 
continuous-time random walk formalism developed by Shlesinger and Landman. 29 The diffusion 
constant for the dimer is computed numerically and plotted in Fig. [71 while the diffusion constant 
for the trimer is given by 

r F 3 H 3 r H 3 F 3 (r F 3 H 3 r H 3 F 3 + r F 3 H 3 r H 3 F 3 ) 

* non AA rot AA non A rot rot rot 1A rot non 1 rat non 1A rot rot 

I r F 3 H 3 r H 3 F 3 r H 3 F 3 + r F 3 H 3 (r F 3 H 3 r H 3 F 3 + (r F 3 H 3 + r H 3 F s )r H 3 F 3 ) 

rot ^ot A1 non A rot rot * non A non A1 rot rot rot 1A non A rot 1 rot non 1A non A rot rot rot 

In Fig.|7J the analytical rates from above have been plotted as diffusion coefficients against tem- 
perature along with experimental data from [12] for the monomer and dimer. The experimental 
barriers, 65+9meV and 73meV for monomer and dimer are both lower than our calculated values, 
though no error bar is given for the dimer experimental barrier. The rate limiting barriers as T— >0K 
are calculated using data at T<20K. The rate limiting barriers, 93meV, 88meV, and 289meV, cor- 
respond to the rate-limiting transition barriers identified above for the monomer, dimer, and trimer. 
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The dimer diffusion slope decreases with increasing temperature due to the influence of both FH 
intermediate states. Higher temperatures samples the FHi ong pathway; this decreases dimer diffu- 
sion as transitions through the FH long state lead to rotation, i.e., FF— >FH long — >HH— >FH long — >FF. 
The prefactors for transitions out of FH short are the lowest for all dimer transitions, and become 
rate limiting at high temperatures. 
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FIG. 7: Experimental results for monomer and dimer compared with analytical diffusion calculations for 
monomer, dimer, trimer, and DFT corrected trimer at low (top) and high (bottom) temperatures. The Arrhe- 
nius fit, in the T— >0K limit, for the monomer, dimer, and trimer systems reflects the rate-limiting diffusion 
barriers 93meV, 88meV, and 289meV respectively. The DFT corrected trimer is calculated by adjusting the 
trimer diffusion barriers to match DFT energy differences and keeping the same prefactors. The experimen- 
tal monomer barrier is 65±9meV from data in the temperature range 19K-25K. The experimental dimer 
barrier is 73meV from data at 24K assuming a prefactor of ITHz. 
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V. DISCUSSION 



The Cu monomer is the basic unit for Cu islands on Ag(l 11), and the correct extrapolation of 
monomer energies and barriers to dimers and trimers indicates the optimized potential is consistent 
with DFT. For the monomer, the site energy difference between FCC and HCP is 14meV, twice 
the difference is seen between the homogeneous dimers FF and HH at 26meV, and two to three 
times the difference in the rot and non trimer pairs at 33meV and 45meV. The diffusion barriers for 
the trimer are also three times that of the monomer 289meV versus 93meV and 247meV versus 
79meV. This linear relationship is explained by the fact that in the trimer diffusional transitions, 
all three atoms move simultaneously over each of their respective bridging sites, thus the trimer as 
a whole experiences a barrier three times as large. In the dimer system, diffusion moves one atom 
at a time and the barrier is comparable to that of the monomer. 

EAM produces higher diffusion barriers for monomer and dimer than in experiment, 1 - but 
gives diffusion barriers that match DFT. This effect is consistent with the observed overestimation 
of surface adsorption energy by DFT calculations. 28 Compared to experiments for monomer and 
dimer, the barriers are overestimated by approximately 10 to 15meV. Since diffusion for both 
the monomer and dimer proceeds one Cu atom at a time, we expect the bridging site between F 
and H to be overestimated by 15meV. For general diffusion barriers, a 15meV reduction should 
be applied for each concurrent Cu atom in the transition when comparing to experiment. For 
example, a three-fold reduction of 45meV will need to be applied to the trimer diffusion barriers. 

Trimers were not included in our optimization database and calculations in DFT and the op- 
timized EAM differ when looking at the energetics between rot and non trimers. DFT calculates 
that the two non trimers are 25meV lower in energy than predicted by our optimized EAM, making 
the ground state trimer configuration F 3 on rather than F 3 ot . We expect the deviation to be mainly 
caused by the center Ag atom under the rot trimer, whose embedded "electron density" is 16% 
higher than an Ag atom in the bulk. This density is beyond the range present in the monomer and 
dimer database. Modifying the embedding function, as done in surface embedded-atom method 
(SEAM) 30 ' 31 may offset this effect by penalizing densities away from the bulk value. Although the 
relative energy between rot and non trimers are not correct, the optimized EAM correctly predicts 
the energy difference between F 3 ot and H 3 ot , and F 3 on and H 3 on . Adding a Cu atom to a rot trimer 
will create a non trimer subsection, and in larger islands, this pairing of rot and non trimers allows 
the correct energy differences to be calculated. We expect the trimer diffusion barrier to remain 
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three times that of the monomer even with the change in ground state. A new estimate of trimer 
diffusion can be computed by splitting the 25meV energy difference between forward and reverse 
diffusion barriers, e.g. lower the rot to non barrier by 12.5meV and raising the non to rot by 
12.5meV. This change does not affect the transition paths and therefore does not change the over- 
all diffusional dynamics of the trimer system, increasing the rate limiting barrier to 292meV from 
289meV (c.f. Fig. [7]). Applying the 45meV over-adsorption correction gives a barrier of 247meV 
for trimer diffusion to compare with experiment. 

VI. CONCLUSION 

We present a method to optimize an EAM potential for heterogeneous surface system using 
ab initio data. The optimized EAM potential reproduces DFT energies for Cu monomers, dimers, 
and trimers on Ag(lll). Diffusion barriers for monomers, dimers, and trimers are calculated to 
be 93meV, 88meV, and 289meV, which match available DFT data, but exceed experimental val- 
ues. To correct for the overestimated barriers, a 15meV reduction is applied for each concurrently 
transitioning Cu atom. We found a 25meV energy discrepancy between rot and non trimers when 
compared with DFT. This discrepancy is not worse for larger islands due to correct energy differ- 
ence between F-trimers and H-trimers calculated by the potential compared with DFT. We expect 
the new EAM potential to accurately describe the diffusion and energetics of larger Cu islands on 
Ag(lll). 

Acknowledgments 

The authors thank John Weaver and Andrew Signor for helpful discussions. This research was 
supported by NSF/DMR grant 0703995, and 3M's Untenured Faculty Research Award. 



1 J.-M. Wen, S.-L. Chang, J. W. Burnett, J. W. Evans, and P. A. Thiel, Phys. Rev. Lett. 73, 2591 (1994). 

2 G. L. Kellogg and A. F. Voter, Phys. Rev. Lett. 67, 622 (1991). 

3 M. C. Bartelt, C. R. Stoldt, C. J. Jenks, P. A. Thiel, and J. W. Evans, Phys. Rev. B 59, 3125 (1999). 

4 G. Antczak and G. Ehrlich, Surface Science Reports 62, 39 (2007). 



14 



5 N. I. Papanicolaou, G. A. Evangelakis, and G. C Kallinteris, Computational Materials Science 10, 105 
(1998). 

6 F. Montalenti and R. Ferrando, Phys. Rev. B 59, 5881 (1999). 

7 A. Bogicevic, Phys. Rev. Lett. 82, 5301 (1999). 

8 H. T. Lorensen, J. K. N0rskov, and K. W. Jacobsen, Phys. Rev. B 60, R5149 (1999). 

9 V. Chirita, E. Munger, J. Greene, and J.-E. Sundgren, Surface Science Letters 436, L641 (1999). 

10 V. Papathanakos and G. A. Evangelakis, Surface Science 499, 229 (2002). 

11 H. Brune, G. S. Bales, J. Jacobsen, C. Boragno, and K. Kern, Phys. Rev. B 60, 5991 (1999). 

12 K. Morgenstern, K.-F. Braun, and K.-H. Rieder, Phys. Rev. Lett. 93, 056102 (2004). 

13 V. Ozolins, C. Wolverton, and A. Zunger, Phys. Rev. B 57, 4816 (1998). 

14 J. C. Hamilton, Phys. Rev. Lett. 77, 885 (1996). 

15 S. M. Foiles, M. I. Baskes, and M. S. Daw, Phys. Rev. B 33, 7983 (1986). 

16 G. Kresse and J. Hafner, Phys. Rev. B 47, RC558 (1993). 

17 G. Kresse and J. Furthmuller, Phys. Rev. B 54, 11169 (1996). 

18 D. Vanderbilt, Phys. Rev. B 41, 7892 (1990). 

19 G. Kresse and J. Hafner, J. Phys.: Cond. Mat. 6, 8245 (1994). 

20 J. P. Perdew and A. Zunger, Phys. Rev. B 23, 5048 (1981). 

21 S. J. Plimpton, J. Comp. Phys. 117, 1 (1995), URL http : //lammps . sandia . gov/index .html 



22 G. Mills and H. Jonsson, Phys. Rev. Lett. 72, 1124 (1994). 

23 G. Henkelman and H. Jonsson, J. Chem. Phys. Ill, 7010 (1999). 

24 G. H. Vineyard, Journal of Physics and Chemistry of Solids 3, 121 (1957). 

25 Y. Mishin, M. J. Mehl, D. A. Papaconstantopoulos, A. F. Voter, and J. D. Kress, Phys. Rev. B 63, 224106 
(2001). 

26 P. L. Williams, Y. Mishin, and J. C. Hamilton, Modelling Simul. Mater. Sci. Eng. 14, 817 (2006). 

27 A. F. Voter, Tech. Rep., Los Alamos Unclassified Technical Report (1993). 

28 C. Stampfl, Catalysis Today 105, 17 (2005). 

29 M. F. Shlesinger and U. Landman, Applied Stochastic Processes (Academic Press, New York, 1980). 

30 M. I. Haftel, Phys. Rev. B 48, 2611 (1993). 

31 M. I. Haftel and M. Rosen, Phys. Rev. B 51, 4426 (1995). 

32 X. Yuan, K. Takahashi, Y. Ouyang, and T Onzawa, J. Phys.: Condens. Matter 15, 8917 (2003). 



15 



